These blocks of code sets base for the subsequent analyses:
- loads all necessary packages
- set working directory
- read in diversity abundance, diversity length, ysi and habitat
- create wide data for Vegan package
- recreate long data with 0 observations
- clean ysi data(remove duplicates)
- create wide temperature data from ysi data
- merge temperature and habitat data with diversity abundance
- assign levels and names to month objects
# load packages
rm(list=ls())
library("tidyverse", lib.loc="~/R/win-library/3.4")
library("vegan", lib.loc="~/R/win-library/3.4")
library("BiodiversityR", lib.loc="~/R/win-library/3.4")
library("knitr", lib.loc="~/R/win-library/3.4")
library("yaml", lib.loc="~/R/win-library/3.4")
library("rmarkdown", lib.loc="~/R/win-library/3.4")
library("magrittr", lib.loc="~/R/win-library/3.4")
library("lubridate", lib.loc="~/R/win-library/3.4")
library("TropFishR", lib.loc="~/R/win-library/3.4")
# make fabsdivlength name of and create working file
#setwd("C:/Users/FABS/Desktop/Rwd_Ben/FABS/working_files/data")
# read in diversity abundnce, diversity length, ysi and habitat data
da <- read.csv ("C:/Users/FABS/Desktop/Rwd_Ben/FABS/working_files/data/da.csv")
dl <- read.csv ("C:/Users/FABS/Desktop/Rwd_Ben/FABS/working_files/data/dl.csv")
habitat <- read.csv ("C:/Users/FABS/Desktop/Rwd_Ben/FABS/working_files/data/habitat.csv")
ysi <- read.csv ("C:/Users/FABS/Desktop/Rwd_Ben/FABS/working_files/data/ysi.csv")
#create da.wide
da.wide <- da %>%
group_by(year, month, day, site, species) %>%
summarise(abundance1 = sum
(abundance, na.rm = TRUE)
)%>%
spread("species", "abundance1")
#replace na with 0
da.wide[is.na(da.wide)] <- 0
#create long data with null observations
da.long <- gather(da.wide, species, abundance, argo:yero)
#remove replicates and duplicates from ysi
ysi <- ysi[!(ysi$replicate == 2), ]
ysi$replicate <- NULL
ysi <- unique(ysi[ ,1:10])
#create temperature wide data from ysi data; can do with salinity and ph also
temp.wide <- ysi
#remove other variables(ph, salinity)
temp.wide[8:10] <- NULL
#remove depth
temp.wide[6] <- NULL
#spread to wide format
temp.wide %<>%
group_by(site, year, month, day, location) %>%
spread("location", "temp")
#merge temperature data
da.wide.habitat.env <- merge(da.wide, temp.wide,
by = c("year", "month", "day", "site"),
all = TRUE)
#merge habitat data
da.wide.habitat.env <- merge(da.wide.habitat.env, habitat, by = "site")
da.wide.habitat.env[5] <- NULL
#assign levels to month
da.wide.habitat.env$month <- factor(da.wide.habitat.env$month,
levels = c(1,2,3,4,5,6,7,8,9,10,11,12),
labels = c("Jan",
"Feb",
"Mar",
"Apr",
"May",
"Jun",
"Jul",
"Aug",
"Sep",
"Oct",
"Nov*",
"Dec"))
da.long$month <- factor(da.long$month,
levels = c(1,2,3,4,5,6,7,8,9,10,11,12),
labels = c("Jan",
"Feb",
"Mar",
"Apr",
"May",
"Jun",
"Jul",
"Aug",
"Sep",
"Oct",
"Nov*",
"Dec"))
#create date variable for TropFishR
#multiple monthly sample events justify bimonthly sample units
dl.month.dup <- dl
dl.month.dup$bimonth <- if_else(dl.month.dup$day < 15, 7, 22)
dl.month.dup$date <- ymd(paste(dl.month.dup$year, dl.month.dup$month, dl.month.dup$bimonth))
#assign levels to month
dl$month <- factor(dl$month,
levels = c(1,2,3,4,5,6,7,8,9,10,11,12),
labels = c("Jan","Feb","Mar","Apr","May","Jun","Jul","Aug","Sep","Oct","Nov","Dec"))
#average lengths
dl.year.month.core <- dl
dl.year.month.core <-
filter(dl.year.month.core, site == "ssp" | site == "hdo" | site == "pba" | site == "chp" | site == "ppo" | site == "wfb") %>%
group_by(year, month, species) %>%
summarise(se = sd(length)/sqrt((length(month))), length = mean(length))
The following presents seasonal abundance and length data for several of the more common species. Abundance for individual species is calculated as the mean of total catch for a site visit(seine sets are lumped together). Length is fork length(FL). Error bars represent standard error of the mean. Summaries are constrained to core sites only.
Length frequency is presented both as absolute (all years and all months) and seasonally (bimonthly). Seasonal mean length is also presented where applicable (only one cohort present).
#filter core sites; generate SE
da.long.year.month.core <- da.long %>%
filter(site == "ssp" | site == "hdo" | site == "pba" | site == "chp" | site == "ppo" | site == "wfb") %>%
group_by(year, month, species) %>%
summarise(se = sd(abundance)/sqrt((length(species))), abundance = mean(abundance))
#length frequency graphics function by year month for individual species
##function of Species and Length filter(outliers)
lfq_year_month <- function(Species, Length) {
plotlfq <- dl.temp <- filter(dl.month.dup, species == Species)
dl.temp <- filter(dl.temp, site == "ssp" | site == "hdo" | site == "pba" | site == "chp" | site == "ppo" | site == "wfb")
dl.temp <- filter(dl.temp, length < Length)
dl.lfq <- lfqCreate(dl.temp, Lname = "length", Dname = "date", Fname = NA, bin_size = 2,
length_unit = "cm", plus_group = FALSE, aggregate_dates = FALSE,
plot = FALSE)
plot(dl.lfq, Fname = "catch", ylab = "Length (mm)")
}
#seasonal mean length function for any species with standard error
#only useful for species with single cohorts; doesn't work for pipefish, perch, herring, etc
mean_length_seasonal <- function(Species) {
filter(dl.year.month.core, species == Species) %>%
ggplot() +
aes(x = month, y = length, colour = factor(year), group = year) +
geom_point(size = 2.5, position = position_dodge(width = 0.25)) +
geom_line(position = position_dodge(width = 0.25)) +
geom_errorbar(aes(ymin = length - se,
ymax = length + se),
width = 0,
size = 1,
position = position_dodge(width = 0.25)) +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
#legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("Mean Length"), colour = "year")
}
#the following is for less frequent species(indicator species)
abundance_seasonal <- function(Species){
filter(da.long.year.month.core, species == Species) %>%
filter(year != "2018") %>%
filter(abundance > 0) %>%
ggplot() +
aes(x = month, y = abundance, colour = species, group = species) +
geom_point(size = 2.5, position = position_dodge(width = 0.25)) +
geom_line( position = position_dodge(width = 0.25)) +
geom_errorbar(aes(ymin = abundance - se,
ymax = abundance + se),
width = 0,
size = 1,
position = position_dodge(width = 0.25))+
facet_wrap(~year)+
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression(Abundance ~ "(mean count/site visit)"))
}
filter(da.long.year.month.core, species == "byro" | species == "cqbr" | species == "unro" | species == "coro" | species == "blro"
| species == "caro" | species == "boro" | species == "qbro" | species == "vero") %>%
filter(year != "2018") %>%
filter(abundance > 0) %>%
ggplot() +
aes(x = month, y = abundance, colour = factor(species, labels = c("Black",
"Black/Yellow",
"Boccacio",
"Canary",
"Copper",
"Copper/Quillback",
"Quillback",
"Unknown",
"Vermillion")), group = species) +
geom_point(size = 2.5, position = position_dodge(width = 0.25)) +
geom_line( position = position_dodge(width = 0.25)) +
geom_errorbar(aes(ymin = abundance - se,
ymax = abundance + se),
width = 0,
size = 1,
position = position_dodge(width = 0.25))+
facet_wrap(~year)+
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
#legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("Rockfish Abundance" ~ "(count/site visit)"), colour = "Species")
The most abundant species are removed to graphically inflate species of lower abundance.
filter(da.long.year.month.core, species == "cqbr" | species == "unro" | species == "coro" | species == "blro"
| species == "caro" | species == "boro" | species == "qbro" | species == "vero") %>%
filter(year != "2018") %>%
filter(abundance > 0) %>%
ggplot() +
aes(x = month, y = abundance, colour = factor(species, labels = c("Black",
"Boccacio",
"Canary",
"Copper",
"Copper/Quillback",
"Quillback",
"Unknown",
"Vermillion")), group = species) +
geom_point(size = 2.5, position = position_dodge(width = 0.25)) +
geom_line( position = position_dodge(width = 0.25)) +
geom_errorbar(aes(ymin = abundance - se,
ymax = abundance + se),
width = 0,
size = 1,
position = position_dodge(width = 0.25))+
facet_wrap(~year)+
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
#legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("Rockfish Abundance" ~ "(count/site visit)"), colour = "Species")
The most abundant species are removed to graphically inflate species of lower abundance.
#rockfish seasonal; byro, cqbr, yellowtail excluded
filter(da.long.year.month.core, species == "unro" | species == "coro" | species == "blro"
| species == "caro" | species == "boro" | species == "qbro") %>%
filter(year != "2018") %>%
filter(abundance > 0) %>%
ggplot() +
aes(x = month, y = abundance, colour = factor(species, labels = c("Black",
"Boccacio",
"Canary",
"Quillback",
"Unknown",
"Vermillion")), group = species) +
geom_point(size = 2.5, position = position_dodge(width = 0.25)) +
geom_line( position = position_dodge(width = 0.25)) +
geom_errorbar(aes(ymin = abundance - se,
ymax = abundance + se),
width = 0,
size = 1,
position = position_dodge(width = 0.25))+
facet_wrap(~year)+
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
#legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("Rockfish Abundance" ~ "(count/site visit)"), colour = "Species")
#length frequency graphics function by year month for individual species
##function of Species and Length filter(outliers)
lfq_year_month("cqbr", 1000)
#execute above function
#only useful for species with single cohorts; doesn't work for pipefish, perch, herring, etc
mean_length_seasonal("cqbr")
#length frequency graphics function by year month for individual species
##function of Species and Length filter(outliers)
lfq_year_month("byro", 1000)
#execute above function
#only useful for species with single cohorts; doesn't work for pipefish, perch, herring, etc
mean_length_seasonal("byro")
These species show no apparent seasonality with respects to numbers or recruitment. There appears to be a similar number from year to year and a similar size distribution.
#from species_individual.R
filter(da.long.year.month.core, species == "enso" | species == "saso" |species == "pasa") %>%
filter(year != "2018") %>%
filter(abundance > 0) %>%
ggplot() +
aes(x = month, y = abundance, colour = factor(species, labels = c("English Sole",
"Sanddab",
"Sand Sole")), group = species) +
geom_point(size = 2.5, position = position_dodge(width = 0.25)) +
geom_line( position = position_dodge(width = 0.25)) +
geom_errorbar(aes(ymin = abundance - se,
ymax = abundance + se),
width = 0,
size = 1,
position = position_dodge(width = 0.25))+
facet_wrap(~year)+
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
#legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("Flatfish Abundance" ~ "(count/site visit)"), colour = "Species")
#length frequency graphics function by year month for individual species
##function of Species and Length filter(outliers)
lfq_year_month("enso", 200)
#execute above function
#only useful for species with single cohorts; doesn't work for pipefish, perch, herring, etc
mean_length_seasonal("enso")
#length frequency graphics function by year month for individual species
##function of Species and Length filter(outliers)
lfq_year_month("saso", 200)
#execute above function
#only useful for species with single cohorts; doesn't work for pipefish, perch, herring, etc
mean_length_seasonal("saso")
#length frequency graphics function by year month for individual species
##function of Species and Length filter(outliers)
lfq_year_month("pasa", 200)
#execute above function
#only useful for species with single cohorts; doesn't work for pipefish, perch, herring, etc
mean_length_seasonal("pasa")
abundance_seasonal("sand")
lfq_year_month("sand", 1000)
abundance_seasonal("ling")
lfq_year_month("ling", 1000)
#execute above function
#only useful for species with single cohorts; doesn't work for pipefish, perch, herring, etc
mean_length_seasonal("ling")
abundance_seasonal("herr")
lfq_year_month("herr", 1000)
abundance_seasonal("shpe")
lfq_year_month("shpe", 1000)
abundance_seasonal("bapi")
lfq_year_month("bapi", 1000)
abundance_seasonal("kegr")
lfq_year_month("kegr", 1000)
abundance_seasonal("tisc")
lfq_year_month("tisc", 100)
abundance_seasonal("busc")
lfq_year_month("busc", 150)
abundance_seasonal("thst")
lfq_year_month("thst", 1000)
abundance_seasonal("stsc")
lfq_year_month("stsc", 1000)
abundance_seasonal("blgo")
lfq_year_month("blgo", 1000)
abundance_seasonal("bago")
lfq_year_month("bago", 1000)
These species show strong annual diffences in abundance and/or recruitment. These could make for potential indicator species. More research into life history traits would be required.
abundance_seasonal("crgu")
lfq_year_month("crgu", 1000)
abundance_seasonal("pegu")
lfq_year_month("pegu", 1000)
abundance_seasonal("grsc")
lfq_year_month("grsc", 1000)
abundance_seasonal("sisc")
lfq_year_month("sisc", 1000)
abundance_seasonal("snpr")
lfq_year_month("snpr", 1000)
#create derived variables for Shannon index and species richness
#new column; shannon H
da.wide.habitat.env$shan <- diversity(
da.wide.habitat.env[(which(colnames(da.wide.habitat.env) == "argo")) :
(which(colnames(da.wide.habitat.env) == "yero"))],
index = "shannon")
#new column; specnumber
da.wide.habitat.env$abundance<- specnumber(
da.wide.habitat.env[(which(colnames(da.wide.habitat.env) == "argo")) :
(which(colnames(da.wide.habitat.env) == "yero"))])
(Kindt’s exact method)
This demonstrates that we have pretty much plateaued in species accumulation. Although, we occasionally find a new species, we only found two new species this year.
specaccum(da.wide.habitat.env[
(which(colnames(da.wide.habitat.env) == "argo")) :
(which(colnames(da.wide.habitat.env) == "yero"))]) %>%
plot(ci.type="polygon", ci.col="lightyellow", xlab = "site visits", ylab = "species number")
The following figures are two metrics (species richness and Shannons index) of diversity filtered, sorted amd displayed in several ways. Richness is calculated as the total number of species observed. Shannonds H take into account both richness and evenness.
Nov*: This sampling event is unique in that it was done late at night and the species comosition is somewhat distinct. Results should be interpreted with caution.
#diversity(spp.richness and shannons) by site and month; all sites
da.month.diversity<- (group_by(da.wide.habitat.env, month, site))%>%
summarise_at(vars(shan:abundance), mean)
#spp.richness
ggplot(da.month.diversity) +
aes(x = month, y = abundance, colour = site, group = site) +
geom_point() +
geom_line() +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1))+
labs(x = "Month", y = expression("Richness" ~ "(species count/site visit)"))
all sites; all years
#shannons
ggplot(da.month.diversity) +
aes(x = month, y = shan, colour = site, group = site) +
geom_point() +
geom_line() +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("Shannon's H" ~ "(Shannons H/site visit)"))
#diversity(spp.richness and shannons) by site and month; core sites
da.month.diversity.core<- (group_by(da.wide.habitat.env, month, site))%>%
summarise_at(vars(shan:abundance), mean) %>%
filter(site == "ssp" | site == "hdo" | site == "pba" | site == "chp" | site == "ppo" | site == "wfb")
core sites; all years
#spp.richness
ggplot(da.month.diversity.core) +
aes(x = month, y = abundance, colour = site, group = site) +
geom_point() +
geom_line() +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
#legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1))+
labs(x = "Month", y = expression("Richness" ~ "(species count/site visit)"))
core sites; all years
#shannons
ggplot(da.month.diversity.core) +
aes(x = month, y = shan, colour = site, group = site) +
geom_point() +
geom_line() +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
#legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("Shannon's H" ~ "(Shannons H/site visit)"))
#diversity(spp.richness and shannons) by year, site and month; all sites
da.year<- (group_by(da.wide.habitat.env, year, month, site))%>%
summarise_at(vars(shan:abundance), mean)
all sites; all years
#spp.richness
ggplot(da.year) +
aes(x = month, y = abundance, colour = site, group = site) +
geom_point() +
geom_line() +
facet_wrap(~year) +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("Richness" ~ "(species count/site visit)"))
all sites; all years
#shannons H
ggplot(da.year) +
aes(x = month, y = shan, colour = site, group = site) +
geom_point() +
geom_line() +
facet_wrap(~year) +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("Shannon's H" ~ "(Shannons H/site visit)"))
#diversity(spp.richness and shannons) by year, site and month; core sites
da.year<- (group_by(da.wide.habitat.env, year, month, site))%>%
summarise_at(vars(shan:abundance), mean) %>%
filter(site == "ssp" | site == "hdo" | site == "pba" | site == "chp" | site == "ppo" | site == "wfb")
core sites; all years
#spp.richness
ggplot(da.year) +
aes(x = month, y = abundance, colour = site, group = site) +
geom_point() +
geom_line() +
facet_wrap(~year) +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
#legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("Richness" ~ "(species count/site visit)"))
core sites; all years
#shannons H
ggplot(da.year) +
aes(x = month, y = shan, colour = site, group = site) +
geom_point() +
geom_line() +
facet_wrap(~year) +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
#legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("Shannon's H" ~ "(Shannons H/site visit)"))
The following are several representations of diverity and abundance by basic habitat type. Habitat was categorized by coarse metrics of substrate and vegetation. Ex) “unal.gravel” = unidentified algae with primarily gravel substrate. Substrate appears to be a good proxy for exposure. Exposed sites have a coarser substrate (sand and gravel), where sheltered sites have a mud (silt and clay) substrate.
Only June, July and August data were used as richness appears to be stable and at a plateau. Richness is calculated as the total number of species observed. Shannonds H take into account both richness and evenness.
summer months (June, July and August)
#diversity by habitat type
#total fish abundance
da.long.habitat <- merge(da.long, habitat, by = "site", incomparables = NA)
group_by(da.long.habitat, year, month, site, veg_subst) %>%
summarise(abundance1 = sum(abundance, na.rm = TRUE)) %>%
group_by(year, month, veg_subst) %>%
summarise(abundance = mean(abundance1), sd(abundance1)) %>%
filter(month == "Aug" | month == "Jul" | month == "Jun") %>%
ggplot() +
aes(x = veg_subst, y = abundance) +
geom_boxplot() +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1))+
labs(x = "vegetation and substrate", y = expression("Fish abundance" ~ "(total fish count/site visit)"))
summer months (June, July and August)
# shannon H by habitat
filter(da.wide.habitat.env, month == "Aug" | month == "Jul" | month == "Jun") %>%
ggplot() +
aes(x = veg_subst, y = shan) +
geom_boxplot() +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1))+
labs(x = "vegetation and substrate", y = expression("Shannons H"))
summer months (June, July and August)
# richness by habitat
filter(da.wide.habitat.env, month == "Aug" | month == "Jul" | month == "Jun") %>%
ggplot() +
aes(x = veg_subst, y = abundance) +
geom_boxplot() +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1))+
labs(x = "vegetation and substrate", y = expression("species count"))
NMDS (Bray-Curtis) ordinations of community composition for the community data.
#extract species observations and sites only
#summarise by site/month (analyses cannot account for day)
#merge habitat data
da.obs <- da.wide.habitat.env
da.obs <- da.obs[1:(which(colnames(da.obs) == "yero"))]
da.obs%<>%
group_by(site, year, month)
da.obs <- summarise_all(da.obs, funs(mean))
da.obs <- merge(da.obs, habitat, by = "site")
NMDS (Bray-Curtis) ordinations of community composition by year in the month of July. Grouped by habitat type. Again, habitat was categorized by coarse metrics of substrate and vegetation. Ex) “unal.gravel” = unidentified algae with primarily gravel substrate. Substrate appears to be a good proxy for exposure. Exposed sites have a coarser substrate (sand and gravel), where sheltered sites have a mud (silt and clay) substrate.
all sites; July
#2014 July NMDS
da.wide.habitat.env.july2014 <- filter(da.obs,
year == "2014",
month == "Jul")
da.wide.habitat.env.july2014 <- column_to_rownames(da.wide.habitat.env.july2014,
var = "site")
nmds.july2014 <- metaMDS(da.wide.habitat.env.july2014[(which(colnames(da.wide.habitat.env.july2014) == "argo")) :
(which(colnames(da.wide.habitat.env.july2014) == "yero"))],
k=3)
ordiplot(nmds.july2014,
display="site",
type="n",
xlim = c(-2, 2))
points(nmds.july2014,
col="black",
pch = (as.integer(da.wide.habitat.env.july2014$veg_subst)))
ordihull(nmds.july2014, da.wide.habitat.env.july2014$veg_subst,
scaling = "symmetric",
lty=(as.integer(da.wide.habitat.env.july2014$veg_subst)))
legend("topright", levels(da.wide.habitat.env.july2014$veg_subst),
pch=1:(length(levels(da.wide.habitat.env.july2014$veg_subst))),
lty=1:(length(levels(da.wide.habitat.env.july2014$veg_subst))))
all sites; July
#2015 July NMDS
da.wide.habitat.env.july2015 <- filter(da.obs,
year == "2015",
month == "Jul")
da.wide.habitat.env.july2015 <- column_to_rownames(da.wide.habitat.env.july2015,
var = "site")
nmds.july2015 <- metaMDS(da.wide.habitat.env.july2015[(which(colnames(da.wide.habitat.env.july2015) == "argo")) :
(which(colnames(da.wide.habitat.env.july2015) == "yero"))],
k=3)
ordiplot(nmds.july2015,
display="site",
type="n",
xlim = c(-2, 2))
points(nmds.july2015,
col="black",
pch = (as.integer(da.wide.habitat.env.july2015$veg_subst)))
ordihull(nmds.july2015, da.wide.habitat.env.july2015$veg_subst,
scaling = "symmetric",
lty=(as.integer(da.wide.habitat.env.july2015$veg_subst)))
legend("topright", levels(da.wide.habitat.env.july2015$veg_subst),
pch=1:(length(levels(da.wide.habitat.env.july2015$veg_subst))),
lty=1:(length(levels(da.wide.habitat.env.july2015$veg_subst))))
all sites; July
#2016 July NMDS
da.wide.habitat.env.july2016 <- filter(da.obs,
year == "2016",
month == "Jul")
da.wide.habitat.env.july2016 <- column_to_rownames(da.wide.habitat.env.july2016,
var = "site")
nmds.july2016 <- metaMDS(da.wide.habitat.env.july2016[(which(colnames(da.wide.habitat.env.july2016) == "argo")) :
(which(colnames(da.wide.habitat.env.july2016) == "yero"))],
k=3)
ordiplot(nmds.july2016,
display="site",
type="n",
xlim = c(-2, 2))
points(nmds.july2016,
col="black",
pch = (as.integer(da.wide.habitat.env.july2016$veg_subst)))
ordihull(nmds.july2016, da.wide.habitat.env.july2016$veg_subst,
scaling = "symmetric",
lty=(as.integer(da.wide.habitat.env.july2016$veg_subst)))
legend("topright", levels(da.wide.habitat.env.july2016$veg_subst),
pch=1:(length(levels(da.wide.habitat.env.july2016$veg_subst))),
lty=1:(length(levels(da.wide.habitat.env.july2016$veg_subst))))
all sites; July
#2017 July NMDS
da.wide.habitat.env.july2017 <- filter(da.obs,
year == "2017",
month == "Jul")
da.wide.habitat.env.july2017 <- column_to_rownames(da.wide.habitat.env.july2017,
var = "site")
nmds.july2017 <- metaMDS(da.wide.habitat.env.july2017[(which(colnames(da.wide.habitat.env.july2017) == "argo")) :
(which(colnames(da.wide.habitat.env.july2017) == "yero"))],
k=3)
ordiplot(nmds.july2017,
display="sites",
type="n",
xlim = c(-2, 2))
points(nmds.july2017,
col="black",
pch = (as.integer(da.wide.habitat.env.july2017$veg_subst)))
ordihull(nmds.july2017, da.wide.habitat.env.july2017$veg_subst,
scaling = "symmetric",
lty=(as.integer(da.wide.habitat.env.july2017$veg_subst)))
legend("topright", levels(da.wide.habitat.env.july2017$veg_subst),
pch=1:(length(levels(da.wide.habitat.env.july2017$veg_subst))),
lty=1:(length(levels(da.wide.habitat.env.july2017$veg_subst))))
NMDS (Bray-Curtis) ordinations of community composition in summer months (June, July and August). Grouped by year and habitat type. Only annual sites were used. Annual sites are those that were visited every year and includes core sites.
annual sites; summer months
#NMDS by year;summer months(June, July, August); all habitats
da.wide.habitat.env.trial <- filter(da.obs,
month == "Jun" | month == "Jul" | month == "Aug",
site == "chp" | site == "fan1" | site == "fan3" | site == "gog1" | site == "hdp" |
site == "kis2" | site == "pba" | site == "ppo" | site == "sni1" | site == "sni2" |
site == "ssp" | site == "wfb")
da.wide.habitat.env.trial$year <- as.factor(da.wide.habitat.env.trial$year)
nmds.trial <- metaMDS(da.wide.habitat.env.trial[(which(colnames(da.wide.habitat.env.trial) == "argo")) :
(which(colnames(da.wide.habitat.env.trial) == "yero"))],
k=3)
ordiplot(nmds.trial,
display="site",
type="n")
points(nmds.trial,
col="black",
pch = (as.integer(da.wide.habitat.env.trial$year)))
ordihull(nmds.trial, da.wide.habitat.env.trial$year,
scaling = "symmetric",
lty=(as.integer(da.wide.habitat.env.trial$year)))
legend("topright", levels(da.wide.habitat.env.trial$year),
pch=1:(length(levels(da.wide.habitat.env.trial$year))),
lty=1:(length(levels(da.wide.habitat.env.trial$year))))
annual sites; summer months
##NMDS by year; summer months(June, July, August); all seagrass
da.wide.habitat.env.trial <- filter(da.obs,
month == "Jun" | month == "Jul" | month == "Aug",
site == "chp" | site == "fan1" | site == "fan3" | site == "gog1" | site == "hdp" |
site == "kis2" | site == "pba" | site == "ppo" | site == "sni1" | site == "sni2" |
site == "ssp" | site == "wfb",
veg_subst == "zostera.mud" | veg_subst == "zostera.sand" | veg_subst == "zostera,gravel")
da.wide.habitat.env.trial$year <- as.factor(da.wide.habitat.env.trial$year)
nmds.trial <- metaMDS(da.wide.habitat.env.trial[(which(colnames(da.wide.habitat.env.trial) == "argo")) :
(which(colnames(da.wide.habitat.env.trial) == "yero"))],
k=3)
ordiplot(nmds.trial,
display="site",
type="n")
points(nmds.trial,
col="black",
pch = (as.integer(da.wide.habitat.env.trial$year)))
ordihull(nmds.trial, da.wide.habitat.env.trial$year,
scaling = "symmetric",
lty=(as.integer(da.wide.habitat.env.trial$year)))
legend("topright", levels(da.wide.habitat.env.trial$year),
pch=1:(length(levels(da.wide.habitat.env.trial$year))),
lty=1:(length(levels(da.wide.habitat.env.trial$year))))
annual sites; summer months
##NMDS by year; summer months(June, July, August); sheltered seagrass
da.wide.habitat.env.trial <- filter(da.obs,
month == "Jun" | month == "Jul" | month == "Aug",
site == "chp" | site == "fan1" | site == "fan3" | site == "gog1" | site == "hdp" |
site == "kis2" | site == "pba" | site == "ppo" | site == "sni1" | site == "sni2" |
site == "ssp" | site == "wfb",
veg_subst == "zostera.mud")
da.wide.habitat.env.trial$year <- as.factor(da.wide.habitat.env.trial$year)
nmds.trial <- metaMDS(da.wide.habitat.env.trial[(which(colnames(da.wide.habitat.env.trial) == "argo")) :
(which(colnames(da.wide.habitat.env.trial) == "yero"))],
k=3)
ordiplot(nmds.trial,
display="site",
type="n")
points(nmds.trial,
col="black",
pch = (as.integer(da.wide.habitat.env.trial$year)))
ordihull(nmds.trial, da.wide.habitat.env.trial$year,
scaling = "symmetric",
lty=(as.integer(da.wide.habitat.env.trial$year)))
legend("topright", levels(da.wide.habitat.env.trial$year),
pch=1:(length(levels(da.wide.habitat.env.trial$year))),
lty=1:(length(levels(da.wide.habitat.env.trial$year))))
annual sites; summer months
##NMDS by year; summer months(June, July, August); exposed seagrass
da.wide.habitat.env.trial <- filter(da.obs,
month == "Jun" | month == "Jul" | month == "Aug",
site == "chp" | site == "fan1" | site == "fan3" | site == "gog1" | site == "hdp" |
site == "kis2" | site == "pba" | site == "ppo" | site == "sni1" | site == "sni2" |
site == "ssp" | site == "wfb",
veg_subst == "zostera.sand" | veg_subst == "zostera.gravel")
da.wide.habitat.env.trial$year <- as.factor(da.wide.habitat.env.trial$year)
nmds.trial <- metaMDS(da.wide.habitat.env.trial[(which(colnames(da.wide.habitat.env.trial) == "argo")) :
(which(colnames(da.wide.habitat.env.trial) == "yero"))],
k=3)
ordiplot(nmds.trial,
display="site",
type="n")
points(nmds.trial,
col="black",
pch = (as.integer(da.wide.habitat.env.trial$year)))
ordihull(nmds.trial, da.wide.habitat.env.trial$year,
scaling = "symmetric",
lty=(as.integer(da.wide.habitat.env.trial$year)))
legend("topright", levels(da.wide.habitat.env.trial$year),
pch=1:(length(levels(da.wide.habitat.env.trial$year))),
lty=1:(length(levels(da.wide.habitat.env.trial$year))))
annual sites; summer months
##NMDS by year; summer months(June, July, August); sandy habitat
da.wide.habitat.env.trial <- filter(da.obs,
month == "Jun" | month == "Jul" | month == "Aug",
site == "chp" | site == "fan1" | site == "fan3" | site == "gog1" | site == "hdp" |
site == "kis2" | site == "pba" | site == "ppo" | site == "sni1" | site == "sni2" |
site == "ssp" | site == "wfb",
veg_subst == "sand")
da.wide.habitat.env.trial$year <- as.factor(da.wide.habitat.env.trial$year)
nmds.trial <- metaMDS(da.wide.habitat.env.trial[(which(colnames(da.wide.habitat.env.trial) == "argo")) :
(which(colnames(da.wide.habitat.env.trial) == "yero"))],
k=3)
ordiplot(nmds.trial,
display="site",
type="n")
points(nmds.trial,
col="black",
pch = (as.integer(da.wide.habitat.env.trial$year)))
ordihull(nmds.trial, da.wide.habitat.env.trial$year,
scaling = "symmetric",
lty=(as.integer(da.wide.habitat.env.trial$year)))
legend("topright", levels(da.wide.habitat.env.trial$year),
pch=1:(length(levels(da.wide.habitat.env.trial$year))),
lty=1:(length(levels(da.wide.habitat.env.trial$year))))
Several summaries of temperature (1 meter depth) filtered, sorted and displayed in several ways.
#create temperature variable
temp.summary<-
group_by(da.wide.habitat.env, year, month, site, veg_subst) %>%
summarise(temp = mean(b2))
error bars are standard error
#monthly temperature by year and habitat type
group_by(temp.summary, year, month, veg_subst) %>%
summarise(se = sd(temp)/sqrt((length(veg_subst))), temp = mean(temp), n()) %>%
filter( year == "2015" | year == "2016" | year == "2017") %>%
ggplot() +
aes(x = month, y = temp, colour = factor(veg_subst, labels = c("macro.sand",
"mud",
"sand",
"ulva.sand",
"unal.cobble",
"unal.gravel",
"unal.mud",
"unal.sand",
"zostera.gravel",
"zostera.mud",
"zostera.sand")), group = veg_subst) +
geom_point(size = 2.5, position = position_dodge(width = 0.25)) +
geom_line( position = position_dodge(width = 0.25)) +
geom_errorbar(aes(ymin = temp - se,
ymax = temp + se),
width = 0,
size = 1,
position = position_dodge(width = 0.25))+
facet_wrap(~year)+
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
#legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1)) +
labs(x = "Month", y = expression("temperature"), colour = "veg_subst")
summer; all sites
#Temperature by habitat type
group_by(temp.summary, year, month, veg_subst)%>%
summarise(se = sd(temp)/sqrt((length(veg_subst))), temp = mean(temp), n()) %>%
filter(month == "Aug" | month == "Jul" | month == "Jun") %>%
ggplot() +
aes(x = veg_subst, y = temp) +
geom_boxplot() +
facet_wrap(~year) +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
legend.position = "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1))+
labs(x = "vegetation and substrate", y = expression("temperature"))
Summer; core sites
error bars are standard error
#temperature by habitat type
filter(temp.summary,
year != "2014",
month == "Jun" | month == "Jul",
site == "ssp" | site == "hdo" | site == "pba" | site == "chp" | site == "ppo" | site == "wfb") %>%
group_by(year, veg_subst)%>%
summarise(se = sd(temp, na.rm = TRUE)/sqrt((length(veg_subst))),
temp = mean(temp, na.rm = TRUE), n()) %>%
ggplot() +
aes(x = veg_subst, y = temp, colour = factor(year, labels = c("2015",
"2016",
"2017")), group = year) +
geom_point(size = 2.5,
position = position_dodge(width = 0.25)) +
geom_errorbar(aes(ymin = temp - se,
ymax = temp + se),
width = 0,
size = 1,
position = position_dodge(width = 0.25)) +
theme_bw() +
theme(axis.text.x = element_text(size = 10, angle = 45, hjust = 1),
#legend.position= "none",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
text = element_text(size = 13, vjust = 0.1))+
labs(x = "vegetation and substrate", y = expression("temperature"), color = "Year")